import numpy as np
from sklearn import datasets, linear_model
import matplotlib.pyplot as plt
from sklearn.metrics import accuracy_score
'''
此文件展示了最初用函数手写实现多层神经网络的程序
其中需要自定义的神经网络参数有：

隐层层数nn.hidden_layer=int(x)
每个隐层的神经元个数nn.hidden_dim=np.array([int(x),int(y),...])
输出层神经元个数nn.output_dim=int(x)
输入层神经元个数是由输入数据决定的，无需自行设置

另外还可以自定义迭代的次数和步长，这两个参数设有默认值:
nn.iteration=2000
nn.step=0.01

在此基础之上我们将程序封装成类的形式，这些将在另一个文件中实现
而使用softmax输出、可视化误分类点、精度对比等要求我们也将在封装成类之后实现
'''
# generate sample data
# 使用的数据集是相互交叉的
#x 是保存坐标的二维array，y 是保存真实标签的一维数组
x, y = datasets.make_moons(200, noise=0.20)
y_true = np.array(y).astype(float) #以float形式保存

# generate nn output target
# matrix t (200,2)
t = np.zeros((x.shape[0], 2))
t[np.where(y==0), 0] = 1 #第一类，y1=1 y2=0
t[np.where(y==1), 1] = 1 #第二类，y1=0 y2=1

class NeuralNetwork():
    iteration=2000 #default iteration number
    step=0.01 #default learning rate
nn=NeuralNetwork()
#set the nn parameters
#根据输入数据x的类型决定输入层的神经元个数
nn.input_dim=np.shape(x)[1]
nn.hidden_layer=3
nn.hidden_dim=np.array([4,4,4])
nn.output_dim=np.shape(x)[1]

nn.W=[]
nn.B=[]
#initialize matrix w and b
#用列表存储输出层和隐层的所有权重和偏置
for layer in range(nn.hidden_layer+1):
    if layer==0:#first hidden layer
        nn.W.append(np.random.randn(nn.input_dim,nn.hidden_dim[0]) \
                / np.sqrt(nn.input_dim))
        nn.B.append(np.zeros((1,nn.hidden_dim[0])))
    elif layer==nn.hidden_layer:#output layer
        nn.W.append(np.random.randn(nn.hidden_dim[layer-1],nn.output_dim) \
                    / np.sqrt(nn.hidden_dim[layer-1]))
        nn.B.append(np.zeros((1,nn.output_dim)))
    else:
        nn.W.append(np.random.randn(nn.hidden_dim[layer - 1], nn.hidden_dim[layer]) \
                      / np.sqrt(nn.hidden_dim[layer - 1]))
        nn.B.append(np.zeros((1, nn.hidden_dim[layer])))

def sigmoid(x):
    return 1.0 / (1 + np.exp(-x))
def forward(nn,x):#输入存储数据的类和输入矩阵
    nn.z=[]
    for layer in range(nn.hidden_layer+1):
        if layer==0:
            nn.z.append(sigmoid(np.dot(x,nn.W[0])+nn.B[0]))
        else:
            nn.z.append(sigmoid(np.dot(nn.z[-1],nn.W[layer])+nn.B[layer]))
    return nn

#show the forward calculation result
forward(nn,x)
y_pred=np.argmax(nn.z[-1],axis=1)
plt.scatter(x[:, 0], x[:, 1], c=y_pred, cmap=plt.cm.Spectral)
plt.title("ini forward")
plt.show()

def BP(nn,x,t):
    for i in range(nn.iteration):
        #正向计算
        forward(nn,x) #在每次前向计算之后会得到这一次的结果nn.z2
        #损失函数和精度计算
        L=0.5*np.sum((nn.z[-1]-t)**2)
        y_pred=np.argmax(nn.z[-1],axis=1)
        acc=accuracy_score(y_true,y_pred)
        if i%100==0:#设置成100次一输出
            print("iteration [%4d] L = %f, acc = %f" % (i, L, acc))

        #反向训练
        delta = []
        for layer in range(nn.hidden_layer,-1,-1):
        #循环的次数应该是hidden_layer+1,可以和前面循环的W B 下标对应
            if layer==nn.hidden_layer: #output layer weight update
                delta.append(nn.z[-1]*(1-nn.z[-1])*(t-nn.z[-1]))
                nn.W[-1]+=nn.step*np.dot(nn.z[nn.hidden_layer-1].T,delta[0])
                nn.B[-1]+=nn.step*np.sum(delta[0],axis=0)
            elif layer==0: #first hidden layer weight update
                delta.insert(0,nn.z[layer] * (1 - nn.z[layer]) * (np.dot(delta[0], nn.W[layer+1].T)))
                nn.W[layer] += nn.step * np.dot(x.T, delta[0])
                nn.B[layer] += nn.step * np.sum(delta[0], axis=0)
            else: #other hidden layers ...
                delta.insert(0,nn.z[layer]*(1-nn.z[layer])*(np.dot(delta[0],nn.W[layer+1].T)))
                nn.W[layer]+=nn.step*np.dot(nn.z[layer-1].T,delta[0])
                nn.B[layer]+=nn.step*np.sum(delta[0],axis=0)

nn.iteration=2000
nn.step=0.01
BP(nn, x, t)
y_pred = np.argmax(nn.z[-1], axis=1)

plt.scatter(x[:, 0], x[:, 1], c=y, cmap=plt.cm.Spectral)
plt.title("ground truth")
plt.show()

plt.scatter(x[:, 0], x[:, 1], c=y_pred, cmap=plt.cm.Spectral)
plt.title("predicted")
plt.show()
